home *** CD-ROM | disk | FTP | other *** search
/ IRIX Base Documentation 2002 November / SGI IRIX Base Documentation 2002 November.iso / usr / share / catman / p_man / cat3 / SCSL / shgeqz.z / shgeqz
Encoding:
Text File  |  2002-10-03  |  11.5 KB  |  265 lines
  1.  
  2.  
  3.  
  4. SSSSHHHHGGGGEEEEQQQQZZZZ((((3333SSSS))))                                                          SSSSHHHHGGGGEEEEQQQQZZZZ((((3333SSSS))))
  5.  
  6.  
  7.  
  8. NNNNAAAAMMMMEEEE
  9.      SHGEQZ - implement a single-/double-shift version of the QZ method for
  10.      finding the generalized eigenvalues  w(j)=(ALPHAR(j) +
  11.      i*ALPHAI(j))/BETAR(j) of the equation  det( A - w(i) B ) = 0  In
  12.      addition, the pair A,B may be reduced to generalized Schur form
  13.  
  14. SSSSYYYYNNNNOOOOPPPPSSSSIIIISSSS
  15.      SUBROUTINE SHGEQZ( JOB, COMPQ, COMPZ, N, ILO, IHI, A, LDA, B, LDB,
  16.                         ALPHAR, ALPHAI, BETA, Q, LDQ, Z, LDZ, WORK, LWORK,
  17.                         INFO )
  18.  
  19.          CHARACTER      COMPQ, COMPZ, JOB
  20.  
  21.          INTEGER        IHI, ILO, INFO, LDA, LDB, LDQ, LDZ, LWORK, N
  22.  
  23.          REAL           A( LDA, * ), ALPHAI( * ), ALPHAR( * ), B( LDB, * ),
  24.                         BETA( * ), Q( LDQ, * ), WORK( * ), Z( LDZ, * )
  25.  
  26. IIIIMMMMPPPPLLLLEEEEMMMMEEEENNNNTTTTAAAATTTTIIIIOOOONNNN
  27.      These routines are part of the SCSL Scientific Library and can be loaded
  28.      using either the -lscs or the -lscs_mp option.  The -lscs_mp option
  29.      directs the linker to use the multi-processor version of the library.
  30.  
  31.      When linking to SCSL with -lscs or -lscs_mp, the default integer size is
  32.      4 bytes (32 bits). Another version of SCSL is available in which integers
  33.      are 8 bytes (64 bits).  This version allows the user access to larger
  34.      memory sizes and helps when porting legacy Cray codes.  It can be loaded
  35.      by using the -lscs_i8 option or the -lscs_i8_mp option. A program may use
  36.      only one of the two versions; 4-byte integer and 8-byte integer library
  37.      calls cannot be mixed.
  38.  
  39. PPPPUUUURRRRPPPPOOOOSSSSEEEE
  40.      SHGEQZ implements a single-/double-shift version of the QZ method for
  41.      finding the generalized eigenvalues w(j)=(ALPHAR(j) +
  42.      i*ALPHAI(j))/BETAR(j) of the equation det( A - w(i) B ) = 0 In addition,
  43.      the pair A,B may be reduced to generalized Schur form: B is upper
  44.      triangular, and A is block upper triangular, where the diagonal blocks
  45.      are either 1-by-1 or 2-by-2, the 2-by-2 blocks having complex generalized
  46.      eigenvalues (see the description of the argument JOB.)
  47.  
  48.      If JOB='S', then the pair (A,B) is simultaneously reduced to Schur form
  49.      by applying one orthogonal tranformation (usually called Q) on the left
  50.      and another (usually called Z) on the right.  The 2-by-2 upper-triangular
  51.      diagonal blocks of B corresponding to 2-by-2 blocks of A will be reduced
  52.      to positive diagonal matrices.  (I.e., if A(j+1,j) is non-zero, then
  53.      B(j+1,j)=B(j,j+1)=0 and B(j,j) and B(j+1,j+1) will be positive.)
  54.  
  55.      If JOB='E', then at each iteration, the same transformations are
  56.      computed, but they are only applied to those parts of A and B which are
  57.      needed to compute ALPHAR, ALPHAI, and BETAR.
  58.  
  59.      If JOB='S' and COMPQ and COMPZ are 'V' or 'I', then the orthogonal
  60.  
  61.  
  62.  
  63.                                                                         PPPPaaaaggggeeee 1111
  64.  
  65.  
  66.  
  67.  
  68.  
  69.  
  70. SSSSHHHHGGGGEEEEQQQQZZZZ((((3333SSSS))))                                                          SSSSHHHHGGGGEEEEQQQQZZZZ((((3333SSSS))))
  71.  
  72.  
  73.  
  74.      transformations used to reduce (A,B) are accumulated into the arrays Q
  75.      and Z s.t.:
  76.  
  77.           Q(in) A(in) Z(in)* = Q(out) A(out) Z(out)*
  78.           Q(in) B(in) Z(in)* = Q(out) B(out) Z(out)*
  79.  
  80.      Ref: C.B. Moler & G.W. Stewart, "An Algorithm for Generalized Matrix
  81.           Eigenvalue Problems", SIAM J. Numer. Anal., 10(1973),
  82.           pp. 241--256.
  83.  
  84.  
  85. AAAARRRRGGGGUUUUMMMMEEEENNNNTTTTSSSS
  86.      JOB     (input) CHARACTER*1
  87.              = 'E': compute only ALPHAR, ALPHAI, and BETA.  A and B will not
  88.              necessarily be put into generalized Schur form.  = 'S': put A and
  89.              B into generalized Schur form, as well as computing ALPHAR,
  90.              ALPHAI, and BETA.
  91.  
  92.      COMPQ   (input) CHARACTER*1
  93.              = 'N': do not modify Q.
  94.              = 'V': multiply the array Q on the right by the transpose of the
  95.              orthogonal tranformation that is applied to the left side of A
  96.              and B to reduce them to Schur form.  = 'I': like COMPQ='V',
  97.              except that Q will be initialized to the identity first.
  98.  
  99.      COMPZ   (input) CHARACTER*1
  100.              = 'N': do not modify Z.
  101.              = 'V': multiply the array Z on the right by the orthogonal
  102.              tranformation that is applied to the right side of A and B to
  103.              reduce them to Schur form.  = 'I': like COMPZ='V', except that Z
  104.              will be initialized to the identity first.
  105.  
  106.      N       (input) INTEGER
  107.              The order of the matrices A, B, Q, and Z.  N >= 0.
  108.  
  109.      ILO     (input) INTEGER
  110.              IHI     (input) INTEGER It is assumed that A is already upper
  111.              triangular in rows and columns 1:ILO-1 and IHI+1:N.  1 <= ILO <=
  112.              IHI <= N, if N > 0; ILO=1 and IHI=0, if N=0.
  113.  
  114.      A       (input/output) REAL array, dimension (LDA, N)
  115.              On entry, the N-by-N upper Hessenberg matrix A.  Elements below
  116.              the subdiagonal must be zero.  If JOB='S', then on exit A and B
  117.              will have been simultaneously reduced to generalized Schur form.
  118.              If JOB='E', then on exit A will have been destroyed.  The
  119.              diagonal blocks will be correct, but the off-diagonal portion
  120.              will be meaningless.
  121.  
  122.      LDA     (input) INTEGER
  123.              The leading dimension of the array A.  LDA >= max( 1, N ).
  124.  
  125.  
  126.  
  127.  
  128.  
  129.                                                                         PPPPaaaaggggeeee 2222
  130.  
  131.  
  132.  
  133.  
  134.  
  135.  
  136. SSSSHHHHGGGGEEEEQQQQZZZZ((((3333SSSS))))                                                          SSSSHHHHGGGGEEEEQQQQZZZZ((((3333SSSS))))
  137.  
  138.  
  139.  
  140.      B       (input/output) REAL array, dimension (LDB, N)
  141.              On entry, the N-by-N upper triangular matrix B.  Elements below
  142.              the diagonal must be zero.  2-by-2 blocks in B corresponding to
  143.              2-by-2 blocks in A will be reduced to positive diagonal form.
  144.              (I.e., if A(j+1,j) is non-zero, then B(j+1,j)=B(j,j+1)=0 and
  145.              B(j,j) and B(j+1,j+1) will be positive.)  If JOB='S', then on
  146.              exit A and B will have been simultaneously reduced to Schur form.
  147.              If JOB='E', then on exit B will have been destroyed.  Elements
  148.              corresponding to diagonal blocks of A will be correct, but the
  149.              off-diagonal portion will be meaningless.
  150.  
  151.      LDB     (input) INTEGER
  152.              The leading dimension of the array B.  LDB >= max( 1, N ).
  153.  
  154.      ALPHAR  (output) REAL array, dimension (N)
  155.              ALPHAR(1:N) will be set to real parts of the diagonal elements of
  156.              A that would result from reducing A and B to Schur form and then
  157.              further reducing them both to triangular form using unitary
  158.              transformations s.t. the diagonal of B was non-negative real.
  159.              Thus, if A(j,j) is in a 1-by-1 block (i.e., A(j+1,j)=A(j,j+1)=0),
  160.              then ALPHAR(j)=A(j,j).  Note that the (real or complex) values
  161.              (ALPHAR(j) + i*ALPHAI(j))/BETA(j), j=1,...,N, are the generalized
  162.              eigenvalues of the matrix pencil A - wB.
  163.  
  164.      ALPHAI  (output) REAL array, dimension (N)
  165.              ALPHAI(1:N) will be set to imaginary parts of the diagonal
  166.              elements of A that would result from reducing A and B to Schur
  167.              form and then further reducing them both to triangular form using
  168.              unitary transformations s.t. the diagonal of B was non-negative
  169.              real.  Thus, if A(j,j) is in a 1-by-1 block (i.e.,
  170.              A(j+1,j)=A(j,j+1)=0), then ALPHAR(j)=0.  Note that the (real or
  171.              complex) values (ALPHAR(j) + i*ALPHAI(j))/BETA(j), j=1,...,N, are
  172.              the generalized eigenvalues of the matrix pencil A - wB.
  173.  
  174.      BETA    (output) REAL array, dimension (N)
  175.              BETA(1:N) will be set to the (real) diagonal elements of B that
  176.              would result from reducing A and B to Schur form and then further
  177.              reducing them both to triangular form using unitary
  178.              transformations s.t. the diagonal of B was non-negative real.
  179.              Thus, if A(j,j) is in a 1-by-1 block (i.e., A(j+1,j)=A(j,j+1)=0),
  180.              then BETA(j)=B(j,j).  Note that the (real or complex) values
  181.              (ALPHAR(j) + i*ALPHAI(j))/BETA(j), j=1,...,N, are the generalized
  182.              eigenvalues of the matrix pencil A - wB.  (Note that BETA(1:N)
  183.              will always be non-negative, and no BETAI is necessary.)
  184.  
  185.      Q       (input/output) REAL array, dimension (LDQ, N)
  186.              If COMPQ='N', then Q will not be referenced.  If COMPQ='V' or
  187.              'I', then the transpose of the orthogonal transformations which
  188.              are applied to A and B on the left will be applied to the array Q
  189.              on the right.
  190.  
  191.  
  192.  
  193.  
  194.  
  195.                                                                         PPPPaaaaggggeeee 3333
  196.  
  197.  
  198.  
  199.  
  200.  
  201.  
  202. SSSSHHHHGGGGEEEEQQQQZZZZ((((3333SSSS))))                                                          SSSSHHHHGGGGEEEEQQQQZZZZ((((3333SSSS))))
  203.  
  204.  
  205.  
  206.      LDQ     (input) INTEGER
  207.              The leading dimension of the array Q.  LDQ >= 1.  If COMPQ='V' or
  208.              'I', then LDQ >= N.
  209.  
  210.      Z       (input/output) REAL array, dimension (LDZ, N)
  211.              If COMPZ='N', then Z will not be referenced.  If COMPZ='V' or
  212.              'I', then the orthogonal transformations which are applied to A
  213.              and B on the right will be applied to the array Z on the right.
  214.  
  215.      LDZ     (input) INTEGER
  216.              The leading dimension of the array Z.  LDZ >= 1.  If COMPZ='V' or
  217.              'I', then LDZ >= N.
  218.  
  219.      WORK    (workspace/output) REAL array, dimension (LWORK)
  220.              On exit, if INFO >= 0, WORK(1) returns the optimal LWORK.
  221.  
  222.      LWORK   (input) INTEGER
  223.              The dimension of the array WORK.  LWORK >= max(1,N).
  224.  
  225.              If LWORK = -1, then a workspace query is assumed; the routine
  226.              only calculates the optimal size of the WORK array, returns this
  227.              value as the first entry of the WORK array, and no error message
  228.              related to LWORK is issued by XERBLA.
  229.  
  230.      INFO    (output) INTEGER
  231.              = 0: successful exit
  232.              < 0: if INFO = -i, the i-th argument had an illegal value
  233.              = 1,...,N: the QZ iteration did not converge.  (A,B) is not in
  234.              Schur form, but ALPHAR(i), ALPHAI(i), and BETA(i), i=INFO+1,...,N
  235.              should be correct.  = N+1,...,2*N: the shift calculation failed.
  236.              (A,B) is not in Schur form, but ALPHAR(i), ALPHAI(i), and
  237.              BETA(i), i=INFO-N+1,...,N should be correct.  > 2*N:     various
  238.              "impossible" errors.
  239.  
  240. FFFFUUUURRRRTTTTHHHHEEEERRRR DDDDEEEETTTTAAAAIIIILLLLSSSS
  241.      Iteration counters:
  242.  
  243.      JITER  -- counts iterations.
  244.      IITER  -- counts iterations run since ILAST was last
  245.                changed.  This is therefore reset only when a 1-by-1 or
  246.                2-by-2 block deflates off the bottom.
  247.  
  248.  
  249. SSSSEEEEEEEE AAAALLLLSSSSOOOO
  250.      INTRO_LAPACK(3S), INTRO_SCSL(3S)
  251.  
  252.      This man page is available only online.
  253.  
  254.  
  255.  
  256.  
  257.  
  258.  
  259.  
  260.  
  261.                                                                         PPPPaaaaggggeeee 4444
  262.  
  263.  
  264.  
  265.